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A linearization technique for optimal 


design of the damping set with 
internal dissipation 


A. Fakharzadeh J.*, H. Alimorad D. and A. Beiranvand 


Abstract 


Considering a damped wave system defined on a two-dimensional domain, 
with a dissipative term localized in an unknown subset with an unknown 
damping parameter, we address the ill-posed shape design problem which 
consists of optimizing the shape of the unknown subset in order to minimize 
the energy of the system at a given time. By using a new approach based 
on the embedding process, first the system is formulated in variational form. 
Then, by transferring the problem into polar coordinates and defining two 
positive Radon measures, we represent the problem in a space of measures. 
Hence, the shape design problem is changed into an infinite linear one whose 
solution is guaranteed. In this stage, by applying two subsequent approxi- 
mation steps, the optimal solution (optimal control, optimal region, optimal 
damping parameter and optimal energy) is identified by a three-phase opti- 
mization search technique. Numerical simulations are also given in order to 
compare this new method with level set algorithm. 


Keywords: Damped wave equation; Dissipation control; Radon measure; 
Search technique; Shape optimization. 
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many technological situations, a given structure whose optimal position 
at rest (for instance), starts to vibrate due to uncontrolled disturbances 
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which we would like to stop. One possibility, although under ideal conditions, 
is described in [16] through damping mechanisms. In the literature, the prob- 
lem of optimal stabilization for the 2-D wave equation has been extensively 
studied from different perspectives (see for instance, [4], [9] and [13]). The 
analysis performed by Hebrard et al. highlights the effect of the over-damping 
phenomenon characteristic of this damped wave equation [2]. Freitas [9] and 
Lopez [13] solved the mentioned problem in which the dissipation vanishes for 
large values of the constant damping coefficient. In 2006, Munch et al. used 
Young measures to solve a similar problem and presented a solution method. 
In that study, the damping coefficient was fixed and the best unknown inter- 
nal region was determined by the use of the gradient descend method [15]. 
In sequence, the best damping coefficient and damping set were determined 
at different times using the level set method [16]. 


In this paper, we solve the problem of finding an optimal observation 
domain w C 2 C R? for general damping wave equations in a new way. We 
optimize not only the placement but also the shape of w, over all possible 
measurable subsets of Q having a certain prescribed measure. Such questions 
are frequently encountered in engineering applications but have rarely been 
treated from the mathematical point of view. In this regard, for the first 
time, we consider a shape optimization problem to find the optimal shape 
and place of a sensor, modeled by a two-dimensional wave equation. The 
objective is to find the shape of the damping set that minimizes the energy 
at some given end time (see [18] and [19]). 


2 Optimal wave damping problem 


Let Q C R? be a domain with piecewise smooth boundary and consider the 
two-dimensional damping wave equation with Dirichlet boundary conditions. 
Consider additionally that w is a subset of Q with positive Lebesgue measure 
which is independent of time t € (0,7). The resulting equation for the 
displacement of the sensor is then ( [2] and [16]) 


Yua = Aya + A(X) Yura = 0, (eb) Ee Nx (0,7), 
Yu,a = 0, (a, t) € 0D x (0,7), (1) 
Yw,a(2, 0) a yo(x), Yu,a(@, 0) = yi (2), DE Q; 


here, a(x) = ax,,(x) € L®©(Q,R*) is a damping function where a € Rt is 
unknown, @ # w is an unknown region in 2, Ow is a smooth and simple closed 
curve boundary which must be identified, x. is a characteristic function of 
w and yo and y; also indicate the initial position and velocity, respectively. 
In addition, regarding the initial conditions, we assume that: 


(yo(x), x(x) € (H?(Q) NH (Q)) x Ho (Q), 
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where H}(Q) (the Sobolev space of order 1 on 2 whose functions are zero on 
the boundary of 2 [22]). System (1) is well-posed [12], whose unique solution 
satisfies: 


Yo,a € C((0,T); H*(Q) 0 Hy (Q) NC*((0,T); Ho (Q) N. C*((0, T); L*(Q)). 


For every t > 0, the energy of system (1) is obtained through ( [16]) 


1 
Bw, 4,t) = 5 f (i alest)P + oa(e,t)Pde, (2) 


which satisfies the following dissipation law (see [23]): 


E(w,a,t) = -| a(x)|Yu,a(a, t)|?dax < 0. 
Q 


Here, Yq denotes the transversal displacement at point x in time t. 

We attempt to find the unknown region w and damping function a(z), si- 
multaneously, through a three-phase optimization procedure which is based 
on an embedding technique. To apply this method, first, we present the 
problem in variational form; next, it is transferred into a new theoretical 
measure problem in which two unknown positive Radon measures in a prod- 
uct space of measures are sought. Then, the solution procedure is explained 
and finally, by a three-phase optimization technique, a nearly optimal shape 
together with the optimal damping function as well as the minimized value 
of system energy are constructed. 

The paper is organized as follows: the next section is devoted to the basic 
deformation in variational form. The aim of Section 4 is to state the problem 
in a polar system. Section 5 deals with the embedding process and approxi- 
mation schemes. In Section 6, based on the previous discussions, we present 
the solution algorithm. Then, two numerical simulations are presented in 
Section 7. Finally, concluding remarks are offered in Section 8. 


3 Basic Deformation 


In general, it is difficult to identify a calssical solution for problem (1); thus 
attempts have usually been made to find a weak (or generalized) solution 
of the problem, which is more applicable in our work. The main idea in 
this replacement is to convert the problem into the variational form. To 
this end, by multiplying the first equation of system (1) with a function 
y € Hj(Q x (0,7)) and using Green’s theorem, to the initial conditions, for 
each t, one obtains: 


J vrcdr— f edyar = | (yet — pas) =0, 
Q Q aq On On 
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[ vdcde =f eayae: 
Q Q 


[ ede f vapar+ | a(x)ypdx = 0. (3) 
Q Q Q 


therefore 


then, we have: 


Integrating both sides of (3) with respect to t over [0,7] gives: 


T T T 
| [ iiedeat - [ | vdedeat + | | e(ahiieanat = 0, (4) 
0 JQ 0 JO 0 Ja 


Double integrating by parts with respect to t from the first left expression 
and integrating the third expression on the left-hand side of (4), we conclude: 


So Ja tiedxdt = foly(P)e(L) - 90) e(0) — yD) 9(T) + y(0)G(0)]dx 


+o So yedadt; . 
Io Sq ow)updadt = fy a(x)[y(T)e(T) - — )— fo yelde 
= fo ala)ly()e(L) — yO)e - Jog ala)ygdadt. 
Now, by substituting the initial conditions of system (1) in (5), we have: 
So Sq tiedadt = fylbP)e(1) — m1 (x)e(0) — y(P) (7) + yo(x)G(0)]dx 
+fo So ybdedt; (6) 
Jo Sq alw)uedadt = fx a(x)[y(T)e(L) — yo(x)e(0)]dax — fy fo alw)yededt. 
By applying (6), equation (4) is changed to: 
fo w(Te(D)de — Ce ie )dx — fy fo yAvdedt + f, a(a)y(T) p(T) dx 
— fo alx)y Oda — [> fo a(x)yodadt + fo yédaxdt (7) 


= foly(z) a = ae 


Moreover, for all (a, t) € OQ x [0, T] by the initial condition, we have y(x,t) = 
0; to apply this condition and using Green’s theorem, we have: 


| y(2, t)y(a, t).ndo = | div(y(a, t)y(a, t))dx = 0. (8) 
aa 


0 


Since the unknown region w must lie in Q and the measure of this unknown 
region must be non-zero, the set of admissible shapes for problem (1) can be 
shown as: 
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Vi = {w C2: |w| = LIQ] }, OL hI (9) 


in which |w| indicates the measure of w and L is a fixed number. This 
constraint can be shown by the following integral relation: 


[acne f ae. (10) 


4 Expressing the problem in polar system 


The mentioned optimal shape design (OSD) problem is defined based on the 
unknown geometrical pair (w, Qw). This pair consists of a measurable set that 
can be regarded as a nonempty region, and a simple closed curve which is its 
boundary. Based on the simplicity property of the curve, our OSD problem 
depends on the geometry which is used. We prefer to solve the appropriate 
problems in polar coordinates since where 0 < @ < 27 and r > O, the curve 
r = r(@) is simple. This simple fact is an essential part in our calculations 
and also in numerical simulations. Hence, let x; = rcos(@) and x2 = rsin(6); 
then, we have: 


— Oyo , 1 ,AY\2 _ Py  1ldy 1 ay 
Vuln = Sp i 72 (59? |, Aylr = Or? r Or | 2 002’ 
and therefore: 
eal Oy.o 1 Oy\5 
Bw, 4,7) = 5 f [hile 8.7) + (GE)? + <a (Gg) Treva, 


Since the nature of has not changed, but rather its representation has 
changed, we use the same symbol and, in the end, the optimal shape is 
shown in polar coordinates. 


Additionally, for every y € Hj(Q x (0,T)), the mentioned constraint in 
(7) can be represented as: 


Jgy¥De(L)rdrdé — fo y(L pas oes - Te Jog yAerdrdédt 
+ fo alr, )y(T)e(L)rdrdé — J. a(r, @)yo(r, Ne i \rdrdo 


~ i a r,) )yprdrdbdt + fo Jag ygrdrdédt = 
(11) 
in which 


= i lyn (r, €)(0) — yo(r, 8) G(0) rar. 
Q 


Moreover, equations (8) and (9) can be represented in polar coordinates as: 
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Jaq Hv(y(r, 4, thy(r, 8, t))rdrd6 = 0, 


12 
ch rdrdé = L ie rdrdé. 2) 


Therefore, the problem of obtaining the optimal shape and damping coeffi- 
cient for minimizing energy in polar coordinates has the following presenta- 
tion: 


Min: E(w,a,T) = 5 Jollg@.8,T)/? + (22)? + 4, (24)? |rdrdé, 

S. to: Jo y(T)y(T)rdrdé — f, y(T)o(L)rdrde — Je Jaq yAyrdrdédt 
+ Jo a(r, A)y(T)p(T)rdrdé — te a(r, 0) yo(r, 8) p(0)rdrdé 
~ fo Jaalr,ygrdrdddt + fy J, yerdrdbdt = &; 
J div(y(r, 9, ty (r, 8, t))rdrdé = 0, Vy € Hy (Q x (0,T)); 
[rdrd0 = L f.,rdrdd. 


(13) 


To solve (12), we change the problem and consider a new one with a different 
formulation. By applying this method, we show how one can obtain the op- 
timal region w, optimal damping function a(x) and the amount of minimized 
energy simultaneously. 


5 Embedding the solution space: metamorphosis 


The solution method which is based on an embedding process involves several 
stages to set up a linear programming problem whose solution converges to 
the solution of the original problem (see [20]). This is one of the outstanding 
advantages of this method even for strongly nonlinear problems. Hence, we 
present a new version of shape measure method to solve the optimal shape 
design (12). First, by defining a new variational formulation, an optimal 
control problem equivalent to the original problem is obtained. Then, a 
measure theoretical approach and a two-stage approximation are used to 
convert the optimal control problem to a finite dimensional LP. The solution 
of this LP is used to construct an approximate solution to the original control 
problem. The proposed approach is practical and accurate enough and its 
accuracy can be improved as far as desired (see [8]). 


5.1 Step 1: Displaying the problem in variational form 


In order to transform the optimal shape design into variational form, we need 
to define some fundamental concepts. The conditions imposed on the func- 
tions and sets will serve two important purposes. First, they are reasonable 
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conditions which are usually met when considering classical problems. Sec- 
ond, they will allow the modification of these classical problems into ones 
which appear to have some advantages over the classical formulation. 


Suppose that J = [0,27] and A; = [0,rq] are the domains of variables 0 
and r,, respectively. For Q = J x Aj, let V C R be a given bounded and 
closed set, and Az = [0,7]. 

Definition 1. Consider the variable r : J° —> Ag as an absolutely continuous 
trajectory function, where J° = (0,27); then, we denote the boundary curve 
of the unknown region w by Ow which is introduced by: 


Ow:r=r(0), dE J. (14) 


Definition 2. By supposing that r is a simple and closed curve and U is a 
bounded closed subset in R, we introduce an artificial control function wu as 
follows: 


u:J° > U, 
u(9) = 7(@) = 9(9,r,u), 


where r(@) € [0,7.,] and the boundary Ow is determined by this variable. 
Definition 3. Let S’ be a bounded closed subset of R and function Y : 
Q x [0,7] > S” be defined in the following way: 


_ Oy 


aS ap 


(15) 
Since y = y(r1, 0) andr; = 1r1(0) € Q, we can write y = y(0) and 6 = 0(r1). In 
this case, we have 0y/Or, = (Oy/00)(d0/dri) and by introducing v : J° > V 
where v = dr) /d0, we can write: 

Oy Oy i 


ar, OO W’ 


where r1(0) € [0,rg] and 00 is determined by this variable. From now on, 
for simiplicity we denote r; also with r by regarding that r; = r € Q; hence 
the derivative of r in w and 2 were shown above by u and v respectively. 


To identify the relationship between r and v as variables, suppose E’ = 
J x A, x V and consider h in C!(E’), then: 


Oh 


Oh 
9g (87) = —(6,r,v)v. (16) 


Or 

In the same way, to display the relationship between functions y and Y as 
two variables, we define E = J x S x S’(where S is the range of y) and 
consider f € C'(E); on the basis of y = y(r,@), it is possible to express 
@ according to y implicitly. Therefore, function f(@,y,Y) can basically be 
displayed as f(0,Y). That is , y is not an independent variable of f; hence, 
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Of /Oy(@,Y) = (Of /00(0,Y)).(06/dy). Now, by regarding (14), we have 
Of /Oy(@,Y) = (Of/00(6,Y)).1/Y) [3]. Since this relationship is satisfied 
for every f € C!(E), it can be concluded that: 


| OF 9 y)rdrad = (8, Y) scrdrad, Vf € C'(E). (17) 


We add this set of constraints to (12) in order to specify the relationship 
between y and Y when they are considered as variables in the problem. 


Since our aim is to identify w and its unknown boundary, we prefer to 
display the constraints of (12) as integrals on the boundary of w as much as 
we can. Thus, the term on the right-hand side of (6), we have: 


A a(x)y(T)o(T) dx — | a(«)4o(«)e(0)de = i: aly(T)y(T) — yo()o(0)) de; 
Qa Q 


Regarding the Green’s theorem : 


N M 


Ow w Ox1 Ore 


suppose M = 0 and N = f5" a(y(T)¢(T) — yo(x)y(0))dx1. Considering the 
fact that a(a) = ax..(x), we have: 


rcos0 
| af ly(T) y(T) — yo(r, 7) p(0)|(*cost — rsinr)(rsint + rcost)drdr. 


(18) 
Also, the area constraint (9) can be presented as: 
1 4 Le 
= r-drd@=L | rdrd§ = L— rdrdé. 
2 Ju re) T Jo Jo 
The independency of the objective function from ¢ results in: 
dL re ees 
E(w,a,T) = or lla r,0,T)|" +(s 2 + a)¥ |rdrdédt. (19) 


Now, by substituting (17) into (10) for every y € H4(Q x (0,T)), we have: 


1 fo fou(T)e(L)rdrdodt — 4 fo fo y(T)o(L)rdrdodt — fy J, yAprardbodt 
+ fo, coe (T)ye(T) — yo(r, 7) e(0)] (Fcost — on + rcost)drdr 
=f JQ ar, *)ygrdrdédt + Lea Jag yerdrdédt = 


Integrating (11), (15) and (16) over [0, 7] implies that: 
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fo Sq Giv(y(r, 8, Heyl, 0, t))rdrdbdt = 0, Ve € H4(Q x (0,T)); 
ape Sq Shrardodt = + L fe fo 26 brdrdodt, Vf € C1(E); 
4 fo So (0,7, v)rdrdodt = 4) f, 22 (6,7, v)urdrdodt, Wh € C1(E’). 


Therefore, problem (12) can be displayed in a new variational form as follows: 


Min: E(w,a,T) = 35 fo Jollg(, 6, T)2 + (& + S)Y?I|rardoat 

S.to: 4 f,,,r2drdd =Ld f” f rdrdodt, 
fo Sq div(y(r, 8, Dy (rr, 4, t))rdrdbdt = 0, Ve € HA(Q x (0,7); 
4) Jo 2£(0,Y )rdrdodt = 4 fy fo 9£$(0,Y )rdrdodt, Vf € C1(E); 
aoe te - r nn =A) f, 26,1, v)urdrdbdt, Wh € C1(E’); 
1S fog(D)y(L)rdrabdt — 4 fo fo y(T)o(L)rdraodt 
ra ee eee ~ fo" fo alr, Oygrdrdbdt; 
+ fon oso” WD) e(T) — yo(0)| (cost — rsinr)(Fsint + rcosr)drdr 


ly Jgygrdrdédt = ®, Vy € Hj(Q x (0,T)). 
(20) 
Now, by determining a suitable control function, we rewrite the problem in 
the form of an optimal control problem. 


5.2 Step 2: Embedding into measure space 


By considering the pair of functions (r,y) as the trajectory and the triple 
(u,v, Y) as the control vector, problem (25) can be considered as an optimal 
control problem. In this manner, we need to present the following definition: 


Definition 4. Quintuplet p = (r,u,v,y,Y) is called admissible when it 
satisfies the following conditions: 


1. The control functions u, vy and Y are bounded and continuous and take 
their values on compact sets U, V and S ; 


2. r is a differentiable function and r(0) = r(27); 
3. y is the bounded solution of the linear damped wave system (1); 
4. The relations (15) and (16) are satisfied. 


The set of all admissible quintuplets is denoted by P. We also suppose that 
P is nonempty; in other words, we suppose that the system is controllable 
(This can be seen in [20], for instance). 
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Let D = [0,T]x Jx Ar xSxS'xV and D’ = Jx A2xU; for any admissible 
quintuplets in P, we define the linear, positive and bounded functionals Ap 
and ['p on C(D) and C(D’) in the following way: 


T(F)= fo fo F(t,6,17,0,y, Y )drdbdt, VF € C(D); 


(21) 
A(G) = Jo,, G(8, 7, u)dd, VG eCc(D); 

Since R® is a locally compact space, by the Heine-Borel theorem ( [21]), 

D C R® is a compact Hausdorff space. Also, for the same reason, D’ is a 

Hausdorff compact space. Therefore, for every given p, Riesz representation 

theorem ( [22]) indicates uniquely two positive Radon measures , .p and Ap, 

so that: 


[p(F) = fp Fdup = up(F), VF €C(D); 


(22) 

Ap(G) = Jy Gd\p = AP(G), VG eEC(D); 
Consequently, any admissible quintuplets can be displayed as (27) by a unique 
pair of measures, say (wp, Ap), in a subset F of M*(D) x M+t(D’), where 
M*(X) is the set of all positive Radon measures on X. Therefore, one can 
transfer problem (25) into a measure space by: 


(r,u,v,y,Y) € P+ > (wp, Ap) € M(D) x M(D’). 


It was proved by Rubio (1986) that such a transformation is an injection. 
To achieve something new, we enlarge the underlying space and consider the 
problem of finding a minimizer pair of measures, say (*, A*), on the space of 
all positive related Radon measures which are just satisfied to the conditions 
of (25) (Not just those that are induced from Riesz Representation theorem). 
Therefore, our method is somehow global. 

We now characterize some properties of admissible pairs. Suppose B is 
an open disc in R? that includes J x Ag; consider C’(B) as the space of 
real-valued continuously differentiable functions on B. Then, for every ¢ in 
C’(B), we define: 


$9(0,7, u) = (0, r)ut+ o6(9,r), V(O,7r,u) € D’. 


Then, since the boundary Ow is a closed and simple curve, we have: 
[ et@.rwdo = f 60,7140 = ds, Yo CB), (23) 
7 J 


where dg = $(27,ra) — 6(0, ra), is still unknown since rq in (0, ra) = (27, ra), 
which is the initial and final point of the closed curve w, is unknown. We will 
explain later that it would be characteristic (see Section (6)). 

Let D(J°) be the space of infinitely differentiable real-valued functions 
with compact support in J°. Define 
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(8,7, u) = r()b(O) + ub), VO,r,u) € D’ 


then, 
[ve r,u)dd = r(2r)w(27) — r(0)¥(0) = 0, Ve Ee D(J°). (24) 


The same situation arises for another special choice of functions in C’(B) 
which are only dependent on variable 0, denoted by Ci(D). Thus, for 
$(0,7r, wu) = v(0), we can have: 


/ V(0)d0 = a,, Vv € Cy(D’), (25) 
J 


where a, is the Lebesgue integral of v over J. 


Regarding the famous properties of admissible quintuplets in P which are 
looked at in (28), (29) and (25), and the definitions of the pair of measures 
(u,A) in (21), problem (25) can now be displayed as follows in which the 
measures A and ys are its unknown variables: 


Min: E(u,d) = wlgzplly(r,, 7)? + (Ge + @)Y?] 
S. to: r(¢9(6,r,u)) =de, Vo € C'(B); 

A958, 7, u)) = 0, Vb € D(J°); 

)=a,, Wy € Ci(D’); 


t))r) =0, sid é RO x (0,T)); 
T)p(T)r) — wey(T)O(L)r) — w(yAgr) — ula(r, O)yor) + w(ysr) 
+A(a Rew) — yoy(0)|(*cost — rsinr)(Fsint + rcost)dr) = ®; 


(26) 
We remind that the theoretical measure problem (26) is linear even though 


the initial problem is highly nonlinear. 


The space M*(D) x M*(D’) is a linear space which will become a locally 
convex topological vector space when it gives the weak*topology. This can be 
defined by the family of semi-norms (ju, A) # |u(F’)| + |A(G)| for F € C(D), 
G € C(D’) and € > 0, which can be on the basis of a family of neighborhoods 
of zero for M+(D) x M*(D’). This family is defined by: 


Uc = {(u,A) € M*(D) x M*(D') : |u(F5)| + A(Gs)| < 69 = 1,2,.... 7} 


12 A. Fakharzadeh J., H. Alimorad D. and A. Beiranvand 


and makes a basis for a weak* topology on the space M*(D) x M+ (D’) (Many 
properties of this topology can be found in the literature such as [24]). In 
this way, M+(D) x M*(D’) under this topology is a Hausdorff space ( [21]). 


The proof of the following theorems can be found in [6], [7] and [20]. 


Theorem 1. 

a) Q C M*(D) x M*(D’) is compact under the weak* topology on M*(D) x 
Mt(D’) 

b) The objective function E(4,A) in problem (26) is continuous. 

c) There exists a pair of measures (u*,*) which are optimal for (26) in the 
set Q C M*t(D) x M*(D’); that is, for every (u,X) € Q, we have: 


E(u*,r") < E(u, A). 


Even though (26) has an optimal solution in Q, it is still very difficult to ob- 
tain the exact solution since the underlying spaces are not finite-dimensional, 
the number of equations is not finite and the unknowns are measures. There- 
fore, it is completely acceptable to seek for a suboptimal solution. Thus, 
first, by choosing suitable dense subsets in the appropriate spaces and then, 
by choosing a finite number of them, the problem is approximated by a semi- 
finite linear programming. 


5.3 Identifying a nearly optimal solution 


It is possible to approximate the solution of problem (26) by the solution of 
a finite-dimensional linear one of sufficiently large dimensions. Besides, by 
increasing the dimension of the problem, the accuracy of the approximation 
can be increased. First, we consider the minimization of (26) not only over 
set Q, but also over its subset called Q(M1, Mo, ..., M7) and defined by only a 
finite number of constraints to be satisfied. This will be achieved by choosing 
countable sets of functions whose linear combinations are dense in appropriate 
spaces and then by selecting a finite number of constraints. Let {¢; :7 € N}, 
{wi :1€ Nb{y:ie€ Nb{g: ie N}, {fi : ie N} and {hi : 2 € 
N} be countable dense (in the topological convergence sense) sets in spaces 
C’(B), D(J°), Ci(D’), H3(Q x (0,T)), C1(E) and C1(E’), respectively. By 
choosing a finite number of functions in each set, the solution of (29) can be 
approximated by the solution of the following one: 
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Min: B(u,d) = w(spllyr,,D)P? + (& + 2)Y%) 


8. 02 O10, 7, 0)) = doy BH 1, 2). 1G; 
Ab? (65750) = 0, FS 1,.2).25,.M5; 
AUs(8)) = 44, 8 = 1,2, 5 Ms: 
A(5r?) — Lu(pr) = 0, én 
p(F;) = 0, ¢=1,2,..., Ma; 
w(G3) — w(H;) = 0, 7 =1,2,...,. Ms; 
w(T;) — Kj) = 0, 7 = 1,2,..., Me; 
(Li 


= 


Of; Ohj 
Hy= p39" = a aeOrvyr, 
Oh; y 
hos apt (G,r, v)vr, L,= eu(T)pi(T)r, 


R= a(r, A)ypir, Ni = yr, 


T; =a 5 ty(T) yi (T) — yo(r, 7) yi (0)|(*cost — rsinr)(*sint + rcost)dr. 
The density property of the selected sets in (27) causes its solution to tend to 
the solution of (26) when Mj, Mo,..., M7 — oo; thus, if numbers Mj, ..., M7 
are selected large enough, (27) is a good approximation of our main problem. 
Now, the number of constraints of the problem is finite, but the problem is 
still infinite since the underlying space is a subspace of measures. It would 
be more convenient if we could approximate the solution just by a solution 
of a simple finite LP. This is precisely our main attention. 


Fakharzadeh et al. (1999) presented that the pair of the optimal measures 
of (25) are in the form of \* = ee Gx, 0(z",) and p* = yng art 0(Z5) 
in which Z* and 2%, belong to dense subsets of D and D’, respectively; 
moreover, 0(t) is a unitary atomic measure with support at the singleton set 
t. Substituting these forms in (27), it might seem that the problem has been 
made even more difficult, since, it is transformed into a non-linear one. But, 
if function E(u, A) can be minimized only with respect to the coefficients a*, 
and (6*,, it will be turned to a linear programming problem. In other words, 
the solution can be obtained approximately by solving just the simple finite 
linear programming like below. If one chooses the points Z* and z>, from 
the dense subsets of D and D’, this fact could be achieved in the second step 
of our approximation. (see [6] for more details): 


14 A. Fakharzadeh J., H. Alimorad D. and A. Beiranvand 


Min: E(a,8,a,ra,) = lau AnO(Zp) 


Bite e ey Oboe = dae: B= Ways 
zy Bills Cn) =O) P28, dig 
Bare BmVs(Zm) = Gs; §=1,2,..., M3; 
Dimer Bm3tm — (PL Dnat Anln = 0; 

Say Olt (Za) 0; i= 1,2,...,Ma: 
ne On GA ay On Z,) = 0, j =1,2,...,Ms; 
wat Onlj(Zn) — WRT! on K;(Zn) = 0, j=1,2,..., Me; 
ey anlLi(Zn) — Pi(Zn) — Qi(Zn) — Ri(Zn) + Ni(Zn)] 
+ EY Bin Ti(Zm) = Bi, Gy ae 
Ayn > 0, MH 1525 NG 
Da on, m = 1,2,...,.M. 

(28) 


Here, we defined 0(Z) = (r/2T)[|y(7,9,T)|? + (1/v? + 1/r?)Y2]. Problem 
(28) is still non-linear because R; and T; are functions of the damping coeffi- 
cient and dy, is unknown since the constant point (0,ra,) of w is unknown. 
Now, by using simultaneous three-phase search techniques for (28), the op- 
timal damping coefficient, r(0) = r(27) = ra,, and the optimal coefficients 
az, +, x, B7,---, 3, would be found as explained in next section. Thus, one 
is able to construct the pair of optimal shape and control function in the 
manner which has been explained in ( [8], [7]). 


6 Algorithm 


To apply the mentioned method for solving problem (1) practically, here we 
present an algorithmic path for the solution procedure. Regarding the pre- 
vious statements, we are able to identify the optimal control and optimal 
region by using the following 4 steps algorithm: 


Step 1: The given sets [0,7], J, Ai,S, S’ and V which form are divided 
into n1,N2,N3,N4,N5 and ng equal parts, and the sets J, Ag and U which 
form w are divided into m,,mz2 and m3 equal parts, respectively; so that, 
the N = nyz.n2.n3.N4.N5.ng, the number of 6-dimensional cells, and the 
M = m,.m2.m3, the number of 3-dimensional cells in the related spaces 
are obtained. Then, in each of these 6-dimensional and 3-dimensional cells 
arbitrary points Z; = (ti,6:,7i, yi, Yi,4i) and z; = (6;,7;,u;) are selected 
respectively. 
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Step 2: For fixed numbers M,, M2, M3, we select M, number of $%(z), 
Mp of w/(z), Mg of v,(z), and for fixed numbers My, M5, Me, Mz, we define 
My, number of F;(Z), Ms of G;(Z) and H;(Z), Mg of I; and K; and M;, of 
Lj, P;, Qi, Ri, N; and T; functions, respectively. Now, one is able to set up the 
finite linear programming (28) with N+ M variables and M, + Mj+...+ M7 
constraints, which is dependent on the variables a and rq. 


Step 3: To solve problem (1), we use an iterative method with two in- 
ner loops and apply a three-phase optimization approach. In the first loop 
by giving a fixing the amount of a, the function Jj, : [0,1] —> R defined 
by Jia(ra) = E*(a, 8,7a) is set up. Then, in the second loop, this function 
is minimized by the use of a standard minimization technique (like a line 
search method) as one of the optimization approaches. We remind that in 
each function, in calculating the standard minimization technique, its related 
LP (28) should be solved (one of the optimization phases). 


If the minimizer of Jiq is called r3, with the optimal value Jf, (ra,) = 
E*(a,8,a,rg,), one is able to set up the function Jz : [0,1] — R by 
Jz(a) = Jj, (rq, ) in the first loop, use a search technique as the last phase of 
the optimization approach, determine the optimal damping coefficient, say 
a*, with the optimal value of energy J3(a*) = E*(a,8,a*,rj,). In this man- 
ner, the damping coefficient, the value of the energy of the system and the 
coefficients a’s and (’s are simultaneously and optimally determined. 


Remark 1. In each stage where alternative optimal cases happen, it suffices 
to select one arbitrarily. 


Step 4: Regarding [20], [8] and [7], for the optimal values aj, a%,...,ay, 
Bi, B3,---, 8% obtained from Step 3, the optimal control and the optimal re- 
gion are determined through the following instructions: 


i) Let 69 = 0 and 6; = 0;-1 + 6} for i = 1,2,...,M. 

ii) For 0 € [0;-1,6;), set u*(0) = u;, where, wu; is the related component. as- 
sociated with point z;. In this manner, according to [20], the nearly optimal 
control can be constructed as a piecewise constant function. 

iii) Let ro = ro, = Ta, using the differential equation u(#) = 7(@), we take 
the following difference equation: 


T; =ri-1+ (0; — O;-1).us, _— 1,2,..,M. 
Therefore, M number points (6;,r;), i = 1,2,...,M of the nearly optimal 


region are determined. Using curve fitting or connecting them by line seg- 
ments, we demonstrate the approximated optimal region. 


Theorem 2. If the used minimization techniques used in Step 3 of the 
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above algorithm are convergent, then the algorithm converges to the optimal 
solution of (1) when M,N, My, Mo,..., M7 tend to infinity. 


Proof. The proof of this theorem is given in Appendix A. 


7 Examples 


Now, to show the efficiency of our method and to explain how it works, we 
solve two numerical examples. It is worth mentioning that these examples 
are taken from [17] and [15] as well as from other studies cited by them in 
order for the readers to be able to compare the two methods. 


Example 1: By defining 2 = (0,1) x (0,1) and selecting a = 10 (con- 
stant), problem (1) was solved by Munch (2009) by the level set method. To 
apply this method, the author used the gradient descend method and also 
supplied some necessity relations by applying the finite difference method. In 
this manner, he used an initial shape to determine the optimal solution. We 
must mention that this approach is very time-consuming and the resulting 
optimal shape is also dependent on the number of iterations as shown on 
Page 25 of [17]. As it is also mentioned in Sub-section 5.1.2 Page 24, the 
results of the problem for variable a, depend on the initial shape. Moreover, 
in this case, the local minima have been obtained which are also completely 
dependent on the initial shape (Page 34 of [17]). 
For the chosen initial conditions: 


yo(9,r) = 100sin(mrcos0)sin(xrsind), 
eee =0, (0r)EQ=IX A), 
The optimal value obtained by Munch (2009) for a = 10 and T = 2 was 
mentioned as E(w,a,T) = 88.17, and for a = 10 and T = 1, as E(w,a,T) = 
249.10. 
We considered the same condition as above, and additional conditions 
that are needed for our method as follows: 


y= l=gP=1)=1 


We supposed L = 0.11164, the area of the unknown region w was equal to 0.7 
and (0,rq,) was a boundary point of w which was determined optimally in 
domain 2 as mentioned in Step 3. Then, by selecting the following functions 
and setting them in (28) for M, = 2, Mo = 10, M3 = 10, M4 = Ms = Mg = 
M7, = 2, we set up the corresponding LP with: 


$9 (0,7, u) = 2rdut+r?; $3 (0,7, u) = 2ru; 
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w(0) = sin(10), 1 =1,2,...,5; Yr (0) = (1 — cos(I'6)), I = 1,2,...,5; 


2n(s—1) 27s Q0r 
Ye fe Sta el GS RE G10. 
“To 10 8 | 10’° 


The functions introduced in (27), which expressed the relationship between 
y and Y and also between r and v, were determined as below: 


F\(0,y,Y) = Oy, Fo(0,y,Y)=@yY, Gi(O,y,Y) =r, 
Go(6,y, 1) = Ory, Ay (8, y, Y) = as A2(8, y,Y) = 2rby, 
hi(0,r,v) = Or, ho(O,r,v) = Orv, 11(0,r,v) =r’, 


In(0,r,v) =r2v, Ky1(0,r,v) =Orv, Ke(0,r,v) = rOv?. 
Also, for i = 1,2 we selected y; = r'sin(i0)t; therefore: 
gi(T) = r'sin(i0), d: = r’sin(i0), 


_ 2’: 1 09: 1 079: 
Agi = + r Or 5 062 


= (i(i —1) +i —7?)r*-*sin(i0)t = 0. 


Thus the functions in (27) were illustrated as: 
L, =r** sin(i0), 
P, =r** sin(i6), 
Q: — 0, 
R; = ayr**1sin(id) 
—_ pithy 4 cos((4% Tcos , cos((i— rcos0 
T; = al! Ute di 2 ( (Gr2)reos#) ((i=2) V4: 


: 4 i+2 1-2 
(rt? 4u7% Pr") 7 sin(i—2)rcos0 sin(i+2)rcos0 
4 ( i—2 i+2 )] 
since: 
T; = a fo y(T)yi(T) — yo(r, 7) yi (0)](*sint + rcost) (cost — rsintr)dr 
oad 


[r’ sin(ir)|(usin T + rcosT)(ucosT — rsint)dr 


=a fo 
=a [ro [rt [r'sin(iT)]|(ru(cos 2T) + (u? — r?)(cos 7) (sin T))dr 
=afy 


cd rity) sin(it) cos 27 + (u2r? — r*+2)(sin(ir) cos T sin 7) ]dr 


cos((4 cos((i—2)Tr ur)*—r'+?) rsin(é— sin(i T1\|rcos 
=a(((= oe (Ci+2)7) 4 cos( 2)7)\) 4 (( = dT (i=2)r (i+2) pee? 


i+2 ; i—2 i—2 i+2 


By dividing each of intervals [0,1] and [0,27] into ten, A; = [0, V2] and 
S' = [—2, 2] into five, S’ = [—2,2] into four, V = [—1,1] and Ag = [0,1] into 
ten and U = [—0.6,0.6] into eleven equal parts, we selected N = 10° points 
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of Z; and M = 1100 points of z;; thus, LP (28) was set up with M+ N 
variables and 30 constraints: 


Min : E(a,8,4,7a,) = Dna On B[L+ (gy + a )¥al 
S104 Sas M Bn (2m Om Urn +42.) = 2r(r4.); 
See Brent = ee : 

eam Bm (Irmcos(l0m) + Umsin(lOm)) = 0, |= 1,2,...,5; 

eM Bm (Urmsin(l’Om) + Um(1 — cos(U/8m))) = 0, U = 1,2,...,5; 
Vinal 3? mOm = 0.7; 

ue An(r2 — OntnUn) = 0; 

Dnal Anan — Onna) = 05 

Dnt Onan — RHE) = 0; 

eo On(TnO2Vn — 2nTnYn) = 0; 

et On (—aynr’ sin(3On)) + ro Bn (bart tim) [S — (eermecete) 


Gan PAE [(sin(1rmCcos0m) — 


+ ike ibis (en(Giaedtrn) 


aa hag —r? sin(20) sin(zrcos@) sin(rrsin6)drdd; 
se a atsin(4By) + EME BB armen U5 


ees a(rS, mtn) (@ [sin(2rmCOSO;n) = Eternal) 


(cos(6rmCosOm ) 
6 


me r? sin(30) sin(zrcos0) sin(rrsinO)drdo; 


Bi + Bo oe = i, 


Boo1 + Bova + -.. + Bi100 = 953 
An >0, n=1,2,...,N; Bm >0, m=1,2,...,M. 
(29) 
For a fixed a = 10, using HBM (Honey-Bee-Method) (see [1]) and the mod- 
ified Simplex method from MATLAB 7.6, we obtained the nearly optimal 
artificial control (see Figure 1), shape (Figure 2), point (0,rj, = 0.35) and 
the also energy value as 175.9. 


Comparing with [17], we found that the optimal value of energy obtained 
by using shape measure method was less, while the obtained optimal regions 
were mostly approximately the same. Additionally, our method took less 
time and the obtained optimal region was independent from the number of 
iterations. 
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Example 2: In this example, we obtained the optimal value of a and w 
simultaneously. As mentioned by Munch et al. (2006), in spite of the ini- 
tial conditions being symmetrical, for large values of a, the obtained optimal 
domain might be non-symmetrical but the value of energy would be less. 

Consider the two-dimensional damped wave equation (1) which is ex- 
pressed on a cyclic domain with center at origin and radius as /2 in time 
interval [0,1]. The aim was to obtain the optimal region w, with a known 
area in the circle such that energy of the system was minimized in the final 
time. The conditions were given same as Example 1, while a was variable 
and supposed to be optimally determined. 

With the same action as Example 1 to set up the related linear program- 
ming (29), the random search method was applied to obtain the optimal value 
of Jig. Also, the optimal value of function Jz was determined by using the 
HBM. In this manner, the obtained results were as follows: 


The optimal damping coefficient a* = 24.3624, rj, = 0.5469 and the value 
of the optimal objective function was 0.101. The optimal (artificial) control 
and the optimal region w are shown in Figures 3 and 4. 

In Section 5.1.4 of [17] for T = 1, a = 29.09 and by 2000 iterations, the 
optimal value of energy was given as 12.56. As emphasized there, for variable 
a, the optimal domain is completely dependent on the initial shape. In our 
method, despite a’s being variable, the obtained optimal region was indepen- 
dent of the initial shape and the amount of optimal energy was considerably 
less, while the time consumed also decreased. 


8 Conclusion 


By doing an embedding process and using the property of positive Radon 
measures, we presented a new and very useful technique for solving the prob- 
lem of minimizing the energy of a damped wave system in an unknown region. 
In this method, the problem was solved by a three-phase optimization search 
technique where the unknown damping coefficient, the region and a point 
of its boundary were found optimally. This method has some advantages 
in comparison to the method used by Munch (2009), since we did not face 
the difficulties mentioned there; such as level set functions being flat, di- 
vergence of the systems with respect to dispersion, and the tendency of time 
toward infinity when damping of numerical waves approaches zero. The most 
important characteristic of our shape measure method is its simplicity and 
its independence from the solution of the initial shape. To obtain the opti- 
mal domain, we just need to use three search techniques while solving linear 
programming problems. Additionally, it is necessary to emphasize that this 
method is much easier, linear and less time-consuming. 
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~0.4b 


~0.6+ 


Figure 1: Optimal control function in Example 1 


Figure 2: Optimal domain with constant damping coefficient a = 10 in Example 1 


0.67 


Figure 3: Optimal control function in Example 2 
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Figure 4: Optimal region in the given domain 2 in Example 2 


Appendix A. Proof of Theorem 2 
To prove this theorem, first, we present the two following lemmas. 


Lemma 1. Consider the linear program (27) consisting of minimizing the 
function up > p(O) = E(u,A) over the set Q(M,,...,M7) of measures in 
M*t(D) x M*(D ) satisfying conditions (27). When M,, Mo,..., M7 tend to 
infinity, 


M,,...,M7)= __ inf Ee 
7(Mui, ..., M7) Q(M,.....Mr) vn 


tends to =infg E(p,A). (This lemma is an extension of Proposition HI.1 
by Rubio (1986)). 


Proof. (i) We prove, first, that the sequence {7(M1, Mo,...,M7)} is con- 
vergent when M,,...,M7 tend to infinity; consider, first, the subsequence 
of (My, Mo,..., M7) in form {n(M,,Mj,...,M,) : My, = 1,2,...}. Since 
Q(M,,...,M1) is a subset of positive Radon measures that satisfies in con- 
straint (26). Therefore, since 

Q(1,1,...,1) D Q(, 2,...,2) D Q(3,3,...,3) D.. D Q(M,..., Mi) D... D Q, 
then, (1, 1,...,1) < (2, 2,...,2) <2. < (M1,..., Mi) <... <9. 


This sequence is non decreasing and bounded above and hence it converges 
to a number ¢ < 7; thus, if « > 0, for M; > N(e), we have: 


In(M,, My, ..., Mi) — ¢| <e (30) 


consider now 7(M,, M2, Mj,...,M1) for both My and Mz larger than N(e). 
Without loss of generality, assume that My > M2. Then: 


(M2, Mo, wy M2) < n(Mi, Mo, Mn, .; M1) < n(M1, My, wy M1); 
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therefore, 
(Mo, Mp, ..., Mz) — ¢ < n(My, Mz, My,..., Mi) — ¢ < n(My, Ma, ..., Mi) — ¢, 
and according to (30), we have: 

In(Mi, Ma,...,.Mi)—¢| <e. (31) 


Now consider (Mi, M2, Ms, M4, wey M1) for M, S Mo = M3 S N(e). Then, 
by the same procedure, one could show that: 


In(M1, M2, M3, Mi, ..., M1) — ¢| Se. (32) 
In a similar manner, for M,; > Mz > M3 > M4, > N(e), we have: 


n(M4, Ma, w+y M4) < (M1, M4, M4, M4, Mi, .-; M1) 
< (My, M2, Ma, Ma, My, wy Mj) < n(M,, Mz, M3, M4,M,, Vs 
< (My, M2, M3, M,, M;,..., M1), 


and by using (30) and (32), we have: 
In( M1, Mz, M3, Ma, M4, ., M1) = ¢| Se. 


Finally, in a similar way, for M,; > Mp > ... > M7 > N(e), one can show 
that: 


n(M7, Mz, w-, M7) < n(M,, Mo, w, M7) < n(M1, Mo, .., Mg, M1), 


and hence: 
In(M1, Mo, .., M7) += ¢| < €. 


Thus, the sequence {7(M1, Mo,...,M7),Mi = 1,2,...,.M7 = 1,2,...} con- 
verges to the number ¢ as Mj,..., M7 tend to infinity. 
(ii) We must prove now that the limit ¢ equals 7 = infg E(w, A). 
We, first, show that this limit ¢ can be computed sequentially. It is known 
that 

C= Pe ee oe La (Mi, ...,Mz)]].--], 
provided that limys,-+00 (M1, ..., Mz) exist since ¢ is a finite number. To 
show the existence of this, we fix M,, Mo,..., Mg and vary M7; since 


Q(M,,..., Me, 1) =) Q(M,,..., Me, 2) 2 caer 2D) Q(M,, ..-; M7) i Mebeed Q; 
thus, 


n(M,, wey Mg, 1) < (My, ..., Me, 2) < er < n(My, ..., M7) < Led < 7). 
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For Mg = 1, 2, ... the non decreasing and bounded above sequence {7(M1, ..., M7), M7 = 
1,2, ...} converges to number ¢(M},..., Mg). Hence, the double limit limyy,_,.5 limyy,+00 
can be computed sequentially. 

Now we define 


Q(M., teey Me) = N7=192(M1, teey Me, M7) = Q(M,, vaky M6, 00). 


For the fixed numbers M1, Mg,..., M5, since ¢(M1,..., Me) = lima, 400 (M1, ..., M7) = 
infgcmy,...,Ms) E(u, A), and 


QM, .., Ms, 1) rnd) Q(M, .-; Ms, Mg) =) Q, 


¢(M,, wey M5, 1) < C(M, wey Ms, 2) < ay < ¢(My,..., Ms, Me) < 1; 
thus ¢(M1,..., Ms, Me) is convergent where Mg — oo and we have: 


C(M,, .-; Ms) = lim C(M, .-; Ms, Me) 
Mbg-co 


(M,, ae Me, M7)). 


= lim [ lim 7 
Meco Mz 00 


In a similar manner, by defining 
Q(M,, .... Ms) = N97,=19(M1, ..., Me), 
we would have: 
¢(M), ..., Ma) = limasyso0 C(M1, ..., Ma, Ms) 
= limp, oo [limas, oo [lim; +00 (M1, ..., Me, M7)]]. 
therefore, in the last stage, we obtain limyy, +. ¢(M1) = ¢. 
(iv) Regarding (i) and (ii), now we can prove ¢ = 7. Let 
Py 4 Oats ad One Mi get M7), 


then P D Q, since Q(M,..., M7) > Q for all My), Mo,..., M7. We can show 
that under the conditions of the problem, Q D P; thus Q = P, that will 
finally imply 
¢= lim ...[ lim ¢(M,...,Me)] =inf E(w, A), 
17 0o Meco Q 
which is the contention in the theorem. 
For this purpose, we prove that if (u,) € P, then they are also in Q. 


For a set of total functions such as ¢,, k = 1,2,..., we have \(¢Z) = dg,, 
according to the definition of P. Based on the definition of total functions, 
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since $9 = ¢-u+¢9 on D’, supp'|br — br, |, SUPp’ |do — ¢o,| tend to zero 
as k tends to infinity, where $f is defined in relationship (27). Therefore: 

|A($9) — dg] = |A($9) — dg — A(OZ) + don 

=| fa,(Gru + b0)d0 — fa (r,U + b0,)48 — (dg — do,)| 

< ee |r = Gr, ||uldO ain fe) lo _ bo, |d0 ae lds _ d¢,| 

< Ky supp’ |or — br,| + Kasupp’|do — $0, | + Kasupp'|o(r, 9) — bx(r, 9), 
tend to zero and hence ($2) = dg, ie. A € Q, where Ky, Ky and K3 are 
constant numbers. Using the similar method, for functions ~9 and vu,, we 
prove the above relationship. Let  € P, thus u(F;) = 0 (F; was defined in 


(27)), based on the definition of total functions, for every given y € H}(D) 
and e€ > 0, there are integer N > 0 and scalars 7;, so that 


N 


suppo,r)supal|div(yyr) — S~ ydiv(yyir)|| < €; 
i=1 


therefore, for every F € C(Q x [0,7]), we have 


M(F)| = |W) — Dil we) 

=|f> fo div(yor)drdodt — TO, fo fo vediv(yyir)drdbat| 

< Ey supio,r)supal|div(ygr) — Swan ydiv(ypir)|| < E.e; (E:constant). 
Because e€ is arbitrary, |j(F')| — 0. Sequentially, by considering the density 
of the functions f; € C(E), h; € C(E ) and vy; € H4(D), we use the same 


method in the case of functions G;, Hj, 1;, Kj, Li, Pi, Qi, Ri, Ni and T; to 
prove that w € Q. Therefore, P C Q and the proof is finished. 


Lemma 2. For every € > 0, the problem of minimizing the function 
aaa Q,O(Z,,) on the set P(M,, M2,...,M7)* described by the inequalities 
(34) has a solution for sufficiently large N = N(e). The solution satisfies: 


N 
n(My, ..., Mz) + ple) < S> anO(Zn) < n(Ma,..., Mr) +6, (33) 


3 
Il 
un 


where p(e) tends to zero as € tends to zero. 


A linearization technique for optimal design of the damping set ... 25 


es re Bm dL (2m) — doy <6, k =1,2,...,; 
=e < OEY Bm! (zm) — 0 < ¢, ay ee 
See =1 1 Bmae(Zm) — Gy SG; s=1,2,..., M3; 
-e< omen 1 a = QE AnTn S€, 
—e< Ory 1 OnFi(Zn) —0<e, 4=1,2,..., Ma; 
BES EEG AF. a, el ee j =1,2,..., Ms; 
=e SONY ond (Zn) — Cnn OnKj(Zn) — 0 < ¢, § = 1,2, 04, Mo3 
-ES ee =1 t On [L i(Zn) = P,(Zn) _ Qi(Zn) = Ri(Zn) ay Ni(Zn)| 

SO See 1 Bm T,(z m) — ®, <€, 7=1,2,..., Mz; 


Aa, >0, n=1,2,..,N; By, >0, m=1,2,...,M. 
This lemma is the developed Theorem III.1 from Rubio (1986). 


Proof. (i) Given € > 0, a set 
{z* k= 1,2,...,My + My + Mg, Z",h =1,2,..., M+ Ms + Mg +My} (35) 


can be introduced, as in the proof of Proposition III.3 of [20], so that inequal- 
ities (34) are satisfied. For sufficiently large N and M, the set 


ON Sheet SD, IN eg = Dg EE 


will contain (37); thus, the set P(M, Mg, ..., M7)* is nonempty for such values 
of N and M, since the N-tuple {6%,k = 1,2,...,Mi + M2 + M3,0,0,..., 0} 
and M-tuple {a7,! = 1,2,..., Ma+...+ M7,0,0,...,0} are in this set. From 
the first set of inequalities of (34), by zm = 1 for all m = 1,2,...,M, we have 
—e< bee 1 8m — At < € and from the fifth set of inequalities for F;(Z,) = 1, 
we have —e€ < ben —1 Qn < ¢€; this set of N and M-tuples with nonnegative 
entries is bounded and also closed; thus, it is compact and the linear function 
psa a,O(Z,,) attains its minimum over this set. Hence, we have: 


N Mi+...+M7 
min S- AnO(Zn) < x a,0(Z*) < n(My,..., M7) +6 (36) 


k=1 
therefore, one of the inequalities of (33) has been proved. For the other, let 
us define Q(M1, ..., Mz)*, by using equations (26) as follow: 
Q(My, .... Mr)® = {ue M*(D), € M*(D)I|NGO{(8,7,4)) — dol S 6 
k=1,2,...,M@i;  |A(H7(0,r,u)) — 0] < e, 1=1,2,..., Ma;..., 
|w(Li) — w(Pi) — w(Qi) — wCRi) + ATi) + WON) — ®i] Se, = 1,2,..., M7}. 
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Then, the set of measures of the type uw = poe An6(Zp) and » = 
eae Bm0O(2Zm) with the coefficients a, and G3, in the set P(M,...,M7)*, is 
a subset of Q(M1,..., M7)*. Thus, 


N 
min > a, O(Zn) = minu(O), (37) 


where the minimum in the left-hand side of this inequality is over the set 
P(M,,...,M7)* and the one in right- hand side is over Q(M1,..., M7)*. Also, 
Q(M,, seen M7) — Nes0Q(Mi, sieiksg M7)‘, and 


QM, iy hg 2 Q(M., 1; Mz) if €1 > €2. (38) 


Let 7(M,...,M7,6) be the infimum of u(©) over the set of measures 
Q(M,,...,M7)*. Then, by (40), we have (M1, ..., M7,61) < n(Mh,..., M7, €2), 
if €, > €. 

It is sufficient for our purposes to consider a sequence of values of « = 1/p 
where p = 1, 2,.... Then, 
1 


1 
(My, ...,Mz,1) < (Mi, ..., Mz, 5) <1. (M,..., M7, re <u. <n, 


the sequence {7(Mj, ..., Mz, ms is non decreasing and bounded above. There- 
fore, it converges to a number y(Mj,..., M7) satisfying 


1 
My,..., Mz) = lim (Mq,....My7,—) = inf @) = (Mj, ..., Mz). 
VM Mr) = eng et gee we) ME 7) 

Thus 
ple) = n(Mi, ..., Mz, €) — (Mi, .., Mz) (39) 


tends to zero as € tends to zero; it follows from (39) and (41) that 
N 
min Sa, (Zx) > minu(®) = n(Mi,..Mz) + o(6): 
n=1 


where the left-hand minimum is over the set P(Mj,...,Mz7)* and the right- 
hand one is over Q(Mj,..., M@7)*. Now we prove Theorem 2 as follow: 


Theorem 2. If the used minimization techniques in Step 3 of the above 
algorithm are convergent, then, the algorithm converges to the optimal solu- 
tion of (1) when M,N, My, Mo,..., M7 tend to infinity. 


Proof. To demonstrate the proof, we have used the proof by contradiction. 
Let (a*, 6*,a*,r*) be the minimizer of E(a,B,a,rq) but (w*,a*) is not the 
minimizer of E(w,a,t), this means that the algorithm does not converge to 
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the solution of (1). Thus, there is (w ,a) such that 
E(w ,a,T) < E(w*,a*,T). (40) 


According to the Riesz representation theorem and considering one to one 
transformation of problem (1), with objective function defined in relationship 
(2), to problem (26), there are unique measures p/’ (@) and p*(®) correspond- 
ing to E(w ,a ,T) and E(w*,a*,T), where: 


# (0) = Bud) <e*() = E(u", 0’). 
According to Lemma 1, we have: 


(My, .... M7) = a inf w(O) > 7 = inf p(@) = mf E(u, A), 
and according to Lemma 2, 
N 
(Mi, a) M7) + ple) < > AnO(Zn) < (Mi, o09 M7) a €, 
therefore, 
w'(0) = inf (0) = inf E(u,d) = BW,X) < EX) 


and 


N 
n +p.) < >> a, 0(Zn) = Bla’, Ba ,7g,) Sn +6, 
n=1 


according to the above relationships: 


thus according to (42), we have: 
E(a’,B arg) < E(a”, B*,a",r9). 


This is in contradiction with what we supposed at the beginning, thus, 
(w*,a*) is minimizer of problem (1). 
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